Defining functions and packages

# 1. Source functions ----

source("functions/graphical_par.R")
source("functions/theme_graph.R")
source("functions/allometry.R")
source("functions/firstup.R")

# 2. Recquired packages ----

library(stringr)
library(plyr)
library(dplyr)
library(tidyr)
library(ggplot2)
library(kableExtra)
library(data.table)
library(formattable)
library(plotly)

# Set theme_perso() as the default ggplot theme
theme_set(theme_graph())

# 3. Load data ----

data_complete <- read.csv("./../data/size_at_age_coral_reef_fishes_data.csv")

Summary

# 1 Summary ----

summary(data_complete)

##            Family               Genus                       Species    
##  Serranidae   :1636   Epinephelus  : 750   Cephalopholis argus  : 518  
##  Acanthuridae :1225   Cephalopholis: 638   Epinephelus merra    : 381  
##  Pomacentridae: 615   Acanthurus   : 489   Myripristis berndti  : 364  
##  Holocentridae: 522   Ctenochaetus : 429   Ctenochaetus striatus: 297  
##  Apogonidae   : 381   Myripristis  : 364   Centropyge flavissima: 251  
##  Balistidae   : 332   Centropyge   : 296   Plectropomus laevis  : 248  
##  (Other)      :1658   (Other)      :3403   (Other)              :4310  
##                   ID            Agei              Ri        
##  GAM18_B099        :  31   Min.   : 0.000   Min.   :0.0000  
##  AC_AC_SC_03_16_001:  28   1st Qu.: 1.000   1st Qu.:0.4031  
##  GAM18_B007        :  28   Median : 4.000   Median :0.6054  
##  GAM18_B123        :  26   Mean   : 4.868   Mean   :0.6959  
##  GAM18_A075        :  25   3rd Qu.: 7.000   3rd Qu.:0.8793  
##  GAM18_A168        :  25   Max.   :30.000   Max.   :3.7845  
##  (Other)           :6206                    NA's   :366     
##      Agecpt            Rcpt             Lcpt             L0p       
##  Min.   : 0.000   Min.   :0.1517   Min.   : 28.11   Min.   :1.450  
##  1st Qu.: 6.000   1st Qu.:0.6445   1st Qu.:145.00   1st Qu.:1.700  
##  Median : 9.000   Median :0.8778   Median :206.13   Median :1.800  
##  Mean   : 9.738   Mean   :1.0064   Mean   :233.01   Mean   :2.021  
##  3rd Qu.:13.000   3rd Qu.:1.1700   3rd Qu.:278.29   3rd Qu.:1.950  
##  Max.   :30.000   Max.   :3.8591   Max.   :984.69   Max.   :4.250  
##                                                                    
##       R0p                Li            Biomass             Location   
##  Min.   :0.01202   Min.   :  1.45   Min.   :    0.4   Gambiers :2004  
##  1st Qu.:0.02780   1st Qu.: 72.39   1st Qu.:   58.0   Hao      : 651  
##  Median :0.03453   Median :150.17   Median :  141.3   Marquesas:1187  
##  Mean   :0.03892   Mean   :160.76   Mean   :  530.1   Moorea   :2112  
##  3rd Qu.:0.04133   3rd Qu.:214.48   3rd Qu.:  313.8   Scilly   : 113  
##  Max.   :0.08564   Max.   :949.58   Max.   :12950.0   Tuamotu  : 302  
##                                     NA's   :603                       
##                                       Observer   
##  Fabien Morat                             :4359  
##  Guillemette de Synéty and Jérémy Wicquart:1368  
##  Jean Bienvenu                            : 642  
##                                                  
##                                                  
##                                                  
## 

# 2. Total number of individuals ----

length(unique(data_complete$ID))

## [1] 848

# 3 Total number of species ----

length(unique(data_complete$Species))

## [1] 51

Individual by species

data_complete %>%
  ddply(., .(Species), summarize, n = length(unique(ID)),
        min_lencap = round(min(Lcpt, na.rm = TRUE), 0),
        max_lencap = round(max(Lcpt, na.rm = TRUE), 0)) %>% 
  kable(., col.names = c("Species", "n", "Min length (TL, mm)", "Max length (TL, mm)")) %>%
  kable_styling(bootstrap_options = c("striped", "hover")) %>%
  column_spec(1, italic = T)

Individuals by location

data_complete %>%
  ddply(., .(Species, Location), summarize, n = length(unique(ID))) %>% 
  spread(., key = Location, value = n) %>% 
  mutate(Gambiers = cell_spec(Gambiers, "html", color = ifelse(is.na(Gambiers), "white", "#446CB3")),
         Hao = cell_spec(Hao, "html", color = ifelse(is.na(Hao), "white", "#446CB3")),
         Marquesas = cell_spec(Marquesas, "html", color = ifelse(is.na(Marquesas), "white", "#446CB3")),
         Moorea = cell_spec(Moorea, "html", color = ifelse(is.na(Moorea), "white", "#446CB3")),
         Scilly = cell_spec(Scilly, "html", color = ifelse(is.na(Scilly), "white", "#446CB3")),
         Tuamotu = cell_spec(Tuamotu, "html", color = ifelse(is.na(Tuamotu), "white", "#446CB3"))) %>% 
  kable(., format = "html", escape = FALSE) %>%
  kable_styling(bootstrap_options = c("striped", "hover")) %>%
  column_spec(1, italic = T)

Agei vs Ri

ggplot(data_complete, aes(x = Agei, y = Ri)) +
  geom_point(color = col_color_graph, fill = col_fill_graph, size = 1, shape = 21)+
  theme(strip.text.x = element_text(size = 8, face = "italic"),
        strip.background = element_rect(colour = "black", fill = col_facet, size = 1),
        plot.title = element_text(colour = col_color_graph),
        plot.subtitle = element_text(colour = "black")) +
  facet_wrap(~Species, scales = "free", ncol = 5) +
  labs(x = "Age (years)", y = "Radius at age i (mm)", title = "Check for outliers", subtitle = "agei vs radi")

Agei vs Li

ggplot(data_complete, aes(x = Agei, y = Li)) +
  geom_point(color = col_color_graph, fill = col_fill_graph, size = 1, shape = 21)+
  theme(strip.text.x = element_text(size = 8, face = "italic"),
        strip.background = element_rect(colour = "black", fill = col_facet, size = 1),
        plot.title = element_text(colour = col_color_graph),
        plot.subtitle = element_text(colour = "black")) +
  facet_wrap(~Species, scales = "free", ncol = 5) +
  labs(x = "Age (years)", y = "Length at age i (mm)", title = "Check for outliers", subtitle = "agei vs li")

Biomass vs Lcpt

# 1. Remove NA from data ----

data_complete_na_rm <- data_complete %>% 
  filter(!is.na(Biomass) & !is.na(Lcpt))

# 2. Get the parameters of the allometric regression for each species ----

coef_allo <- ddply(data_complete_na_rm, .(Species), function(x) allometry(data = x, xreg = "Lcpt", yreg = "Biomass", graph = FALSE))

# 3. Generate predicted values ----

pred_values <- ddply(data_complete_na_rm, .(Species), function(x) seq(0, max(x$Lcpt), length.out = 100)) %>% # Create new x
               gather(., pos, xnew, -1) %>% # Gather multiple columns into two columns
               select(Species, xnew) %>% # Delete useless column
               arrange(Species) %>% # Sort DF by species
               merge(., coef_allo, by = "Species") %>% # Merge the DF with parameters of allometric regression
               mutate(ynew = a*xnew^b) # Calculate new y

# 4. Plot ----

ggplot() +
  geom_line(data = pred_values, aes(x = xnew, y = ynew), color = col_color_graph) +
  geom_point(data = data_complete, aes(x = Lcpt, y = Biomass), color = col_color_graph, fill = col_fill_graph, size = 1, shape = 21)+
  theme(strip.text.x = element_text(size = 8, face = "italic"),
        strip.background = element_rect(colour = "black", fill = col_facet, size = 1),
        plot.title = element_text(colour = col_color_graph),
        plot.subtitle = element_text(colour = "black")) +
  facet_wrap(~Species, scales = "free", ncol = 5) +
  labs(x = "Length at capture (mm)", y = "Biomass at capture (g)", title = "Check for outliers", subtitle = "lencap vs biomass")

Export table for Morat et al

# 1. Get the type of each variable ---- 

table_type <- data_complete %>% 
  summarise_all(class) %>% 
  gather(.) %>% 
  rename(Variable = key, Type = value) %>% 
  mutate(Type = firstup(Type), # First letter in uppercase
         Column = 1:nrow(.)) # Add the numbering of the variable

# 2. Get the range of each numeric variable ----

table_range <- data_complete %>% 
  summarise_if(is.numeric, ~paste(round(min(.x, na.rm = TRUE), digits = 3), 
                                  round(max(.x, na.rm = TRUE), digits = 3), 
                                  sep = " - ")) %>% 
  gather(.) %>% 
  rename(Variable = key, Range = value)

# 4. Merge the three datasets and rename variables ----

merge(table_type, table_range, by = "Variable", all.x = TRUE) %>% 
  merge(., read.csv2("./../data/variables_description_for_Morat_et_al.csv"), by = "Variable", all.x = TRUE) %>% 
  select(Column, Variable, Description, Unit, Type, Range) %>% 
  arrange(Column) %>% 
  write.csv2(., "./../data/variables_complete_for_Morat_et_al.csv", row.names = FALSE)

Number of NA by variable

data_complete %>%
  summarise_all(~(sum(is.na(.)))) %>% 
  t(.) %>% 
  as.data.frame() %>% 
  rename("NA" = "V1")

##           NA
## Family     0
## Genus      0
## Species    0
## ID         0
## Agei       0
## Ri       366
## Agecpt     0
## Rcpt       0
## Lcpt       0
## L0p        0
## R0p        0
## Li         0
## Biomass  603
## Location   0
## Observer   0

Abberant trajectories

# 1. Plotly for Acanthurus lineatus ---- 

data_complete %>% 
  filter(Species == "Acanthurus lineatus") %>% 
  plot_ly(., x = ~Agei, y = ~Li, text = ~paste("ID:", ID))

# 2. Plotly for Lutjanus kasmira ---- 

data_complete %>% 
  filter(Species == "Lutjanus kasmira") %>% 
  plot_ly(., x = ~Agei, y = ~Li, text = ~paste("ID:", ID))

Reproducibility

# 1. Reproducibility ----

sessionInfo()

## R version 3.6.0 (2019-04-26)
## Platform: x86_64-w64-mingw32/x64 (64-bit)
## Running under: Windows 10 x64 (build 17763)
## 
## Matrix products: default
## 
## locale:
## [1] LC_COLLATE=French_France.1252  LC_CTYPE=French_France.1252   
## [3] LC_MONETARY=French_France.1252 LC_NUMERIC=C                  
## [5] LC_TIME=French_France.1252    
## 
## attached base packages:
## [1] stats     graphics  grDevices utils     datasets  methods   base     
## 
## other attached packages:
##  [1] plotly_4.9.0        formattable_0.2.0.1 data.table_1.12.2  
##  [4] kableExtra_1.1.0    ggplot2_3.1.1       tidyr_0.8.3        
##  [7] dplyr_0.8.1         plyr_1.8.4          stringr_1.4.0      
## [10] RColorBrewer_1.1-2  extrafont_0.17     
## 
## loaded via a namespace (and not attached):
##  [1] tidyselect_0.2.5  xfun_0.7          purrr_0.3.2      
##  [4] colorspace_1.4-1  htmltools_0.3.6   viridisLite_0.3.0
##  [7] yaml_2.2.0        rlang_0.3.4       pillar_1.4.1     
## [10] later_0.8.0       glue_1.3.1        withr_2.1.2      
## [13] munsell_0.5.0     gtable_0.3.0      rvest_0.3.4      
## [16] htmlwidgets_1.3   evaluate_0.14     labeling_0.3     
## [19] knitr_1.23        httpuv_1.5.1      crosstalk_1.0.0  
## [22] Rttf2pt1_1.3.7    highr_0.8         Rcpp_1.0.1       
## [25] xtable_1.8-4      readr_1.3.1       promises_1.0.1   
## [28] scales_1.0.0      webshot_0.5.1     jsonlite_1.6     
## [31] mime_0.6          hms_0.4.2         digest_0.6.19    
## [34] stringi_1.4.3     shiny_1.3.2       grid_3.6.0       
## [37] tools_3.6.0       magrittr_1.5      lazyeval_0.2.2   
## [40] tibble_2.1.1      crayon_1.3.4      extrafontdb_1.0  
## [43] pkgconfig_2.0.2   xml2_1.2.0        assertthat_0.2.1 
## [46] rmarkdown_1.13    httr_1.4.0        rstudioapi_0.10  
## [49] R6_2.4.0          compiler_3.6.0

Jeremy WICQUART | jeremywicquart@gmail.com | 2019-07-15 16:43:35